;EFA  
;by Rubén Ledesma
;e-mail rdledesma@gmail.com

;*SLOTS

(defmeth efa-model-object-proto :data
  (&optional (values nil set))
  (if set (setf (slot-value 'data) values))
  (slot-value 'data))

(defmeth efa-model-object-proto :nvars
  (&optional (values nil set))
  (if set (setf (slot-value 'nvars) values))
  (slot-value 'nvars))

(defmeth efa-model-object-proto  :nobs
  (&optional (values nil set))
  (if set (setf (slot-value 'nobs) values))
  (slot-value 'nobs))

(defmeth efa-model-object-proto  :niter
  (&optional (values nil set))
  (if set (setf (slot-value 'niter) values))
  (slot-value 'niter))


(defmeth efa-model-object-proto :pca-or-fa
  (&optional (values nil set))
  (if set (setf (slot-value 'pca-or-fa) values))
  (slot-value 'pca-or-fa))

(defmeth efa-model-object-proto :extrac-communalities
  (&optional (values nil set))
  (if set (setf (slot-value 'extrac-communalities) values))
  (slot-value 'extrac-communalities))

(defmeth efa-model-object-proto :initial-communalities
  (&optional (values nil set))
  (if set (setf (slot-value 'initial-communalities) values))
  (slot-value 'initial-communalities))

(defmeth efa-model-object-proto :init-eigenvalues
  (&optional (values nil set))
  (if set (setf (slot-value 'init-eigenvalues) values))
  (slot-value 'init-eigenvalues))
(defmeth efa-model-object-proto :corcov-value
  (&optional (values nil set))
  (if set (setf (slot-value 'corcov-value) values)); 0=varianzas 1=correlaciones
  (slot-value 'corcov-value))

(defmeth efa-model-object-proto :method-value
  (&optional (values nil set))
  (if set (setf (slot-value 'method-value) values)); 0=simulate  1=permute
  (slot-value 'method-value))

(defmeth efa-model-object-proto :nfactors
  (&optional (values nil set))
  (if set (setf (slot-value 'nfactors) values))
  (slot-value 'nfactors))

(defmeth efa-model-object-proto :factor-load
  (&optional (values nil set))
  (if set (setf (slot-value 'factor-load) values))
  (slot-value 'factor-load))

(defmeth efa-model-object-proto :extrac-eigen
  (&optional (values nil set))
  (if set (setf (slot-value 'extrac-eigen) values))
  (slot-value 'extrac-eigen))

;*OPTIONS* 
(defmeth efa-model-object-proto :options ()

  (setf corcov (send choice-item-proto :new 
                            (list "Covariances"
                                  "Correlations")
                            :value 1))

  (setf method (send choice-item-proto :new 
                    (list                              
                     "Least-squares"
                     "Maximum-likelihood"
                     "Principal-Axis"
                     "Principal-Component")
                            :value 0))

  (setf text-nfactors (send edit-text-item-proto :new "3" :text-length 5))
  (setf text-niter  (send edit-text-item-proto :new "100" :text-length 5))

  (setf OK (send modal-button-proto :new "  Ok"
                 :action
               #'(lambda ()
                   (let (
                         (dialog (send ok :dialog))
                         )
                     (send self :corcov-value (send corcov :value))
                     (send self :method-value (send method :value))
                     (send self :niter (read-from-string (send text-niter :text)))
                     (send self :nfactors (read-from-string (send text-nfactors :text)))   
                     ))))
                     
  (setf cancel (send modal-button-proto :new "  Cancel"
                     :action
               #'(lambda ()
                   (let (
                         (dialog (send cancel :dialog))
                         )
                     (send dialog :modal-dialog-return nil)))))

(setf vista-efa-dialog
      (send modal-dialog-proto :new 
            (list 
            (list "Compute from:       " corcov)
             (list "Extraction method: " method)            
             (list "Number of Factors  " text-nfactors)
             (list "Max iter           " text-niter)
             (list ok cancel))
			:default-button ok :title "Options for EFA"
            ))
(setf result (send vista-efa-dialog :modal-dialog))
             result)

;*ANALISIS

(defmeth efa-model-object-proto :analysis ()   
  (let* (
         (data (send current-data :active-data-matrix '(numeric)))
         (nobs (send current-data :active-nobs))
         (nvars (send current-data :active-nvar '(numeric)))
         (correlation-matrix (correlation-matrix data))
         (init-eigenvalues (eigenvalues correlation-matrix))
         (Extraction (select (list 'least-squares-factor-analysis 'maximum-likelihood-factor-analysis 'paxis-factor-analysis 'pca-factor-analysis)  (send self :method-value)))
         )
    (setf EFA (funcall extraction correlation-matrix (send self :nfactors)))
    (send self :factor-load (first efa))   
    (setf estimated-uniqueness (second efa))
    (send self :init-eigenvalues init-eigenvalues)
    (setf initial-uniqueness  (uniqueness-bound correlation-matrix))
    (send self :initial-communalities (if ( = (send self :method-value) 3) (diagonal correlation-matrix) (- 1 initial-uniqueness)))
    (send self :extrac-communalities (- 1 estimated-uniqueness))
    (send self :data data)
    (send self :nvars nvars)
    )
  )

(defmeth efa-model-object-proto :Report
  (&key (stream t) 
        (dialog nil))
  ;(if (not (eq current-object self)) (setcm self))
  (let* ((w nil)
         (data (send self :data))
         (var-labels (send current-data :variables))
         (nvars  (send self :nvars))
         (nobs (send self :nobs))
         (nfactors (send self :nfactors))
         (factor-load (send self :factor-load))
         (Factor-Labels (mapcar #'(lambda (x) (format nil " Factor ~a  " x))(+ 1 (iseq nvars))))
         (extrac-eigen (map-elements 'sum (column-list  (* factor-load (combine factor-load)))));; suma de las saturaciones al cuadrado para cada factor    
         )    
    (setf w (report-header (send self :title "EFA Report")))           
    (display-string (format nil "~%EXPLORATORY FACTOR ANALYSIS PLUG-IN ~%by Ruben Ledesma~%Extraction method adapted from deLeew (1995) ~2%MODEL: ~a~2%" (send self :name)) w)    
    (display-string (format nil "VARIABLES: ~a~2%" var-labels ) w)
    (display-string (format nil "EXTRACTION METHOD: ~a~2%" (select (list 'least-squares 'maximum-likelihood 'paxis-factor-analysis 'pca-factor-analysis) (send self :method-value))) w)
    
 (display-string  (format nil "~%Correlation matrix:~%") w)
    (print-matrix-to-window (correlation-matrix data)  w  :row-heading "" :row-labels var-labels :column-heading "" :column-labels var-labels  :decimals 3)

   (display-string (format nil "~%Determinant: ~a~2%" (determinant (correlation-matrix data))) w)

    (display-string  (format nil "~%Initial Eigenvalues:~%") w)
    (print-matrix-to-window  (transpose (make-array (list 3 nvars) :initial-contents (combine (send  self :init-eigenvalues) (/ (send self :init-eigenvalues) nvars) (cumsum (/ (send self :init-eigenvalues) nvars)))))  w :row-heading "" :column-heading "" :row-labels factor-Labels :column-labels (list "E-value" "Prop" "CumProp") :decimals 3)   
        
    (display-string (format nil "~%Extracted Eigenvalues:~%") w)
    (print-matrix-to-window  (transpose (make-array (list 3 nfactors) :initial-contents  (combine extrac-eigen  (/ extrac-eigen nvars) (cumsum (/ extrac-eigen nvars)))))  w :row-heading "" :column-heading "" :row-labels (select factor-Labels  (iseq nfactors)) :column-labels (list "E-value" "Prop" "CumProp") :decimals 3)   

    (display-string (format nil "~%Variables loading:~%") w)
    (print-matrix-to-window  factor-load  w :row-heading "" :column-heading "" :row-labels var-labels  :column-labels (select factor-Labels  (iseq nfactors)) :decimals 3)   

    (display-string (format nil "~%Communalities:~%") w)

    (print-matrix-to-window  (transpose (make-array (list 2 nvars) :initial-contents  (combine  (send self :initial-communalities) (send self :extrac-communalities)))) w :row-heading "" :column-heading "" :row-labels var-labels  :column-labels (list "Initial" "Extraction") :decimals 3)

   (display-string (format nil "~%Reproduced correlation matrix:~%") w)

    (print-matrix-to-window   (matmult factor-load (transpose factor-load))  w :row-heading "" :column-heading "" :row-labels var-labels  :column-labels var-labels :decimals 3)

   (display-string (format nil "~%Residual correlation matrix:~%") w)

    (print-matrix-to-window   (- (correlation-matrix data) (matmult factor-load (transpose factor-load)) ) w :row-heading "" :column-heading "" :row-labels var-labels  :column-labels var-labels :decimals 3)


    (send w :fit-window-to-text)
    )
  )


(defun paxis-factor-analysis 
  (c p &key (epsilon 1e-6) (nmax 200) (verbose t))
"Args: (c p)
Exploratory Principal Axis Factor Analysis of the covariance
or correlation matrix C. Extracts p factors. Starts with SMC uniqueness
estimates."
(let* (
       (u (uniqueness-bound c))
       (d (- c (diagonal u)))
       (a (truncated-eigen-decomposition d p))
       (s (matmult a (transpose a)))       
       (it 0)
       )
(loop
(setf it (1+ it))
(setf u (- (diagonal c) (diagonal s)))
(setf d (- c (diagonal u)))
(setf a (truncated-eigen-decomposition d p))
(setf s (matmult a (transpose a)))
(if   (= it nmax)
    (return  (list a u))
    (setf a a))
)))


(defun pca-factor-analysis 
  (c p)
"Args: (c p)
Principal Component Analysis of the covariance
or correlation matrix C. Extracts p components."
(let* (       
       (a (truncated-eigen-decomposition c p))
       (u (- 1 (diagonal (matmult a (transpose a)))))
       )
  (list a u)))

*;Create Data

(defmeth efa-model-object-proto :create-data 
  (&key (dialog nil)
        )
  
  (if (not (eq current-object self)) (setcm self))
  
  (let ((creator (send *desktop* :selected-icon))
        )
    (data (strcat "Score-" (send self :name))
          :created creator
          :creator-object self
          :title (strcat "Factor scores for " (send self :title))
          :data  (combine (normalize (matmult (normalize (send self :data)) (send self :factor-load))))
          :variables (mapcar #'(lambda (x) (format nil "Fact~a" x)) (iseq 1 (send self  :nfactors)))                    :types (repeat "Numeric" (iseq (send self :nfactors))))))

